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I. INTRODUCTION 

The calculation of the frequencies of a drum of arbi- 
trary shape is a problem of formidable difficulty. In fact, 
the exact solutions to the Helmholtz equation arc known 
only in special cases, namely for the equilateral triangle, 
for rectangles and for circles and ellipses, where they are 
expressed in terms of trigonometric, Bessel and Mathieu 
functions. For general shapes, including those obtained 
performing small perturbations of the domains just men- 
tioned, exact results are not available. In these cases 
one is forced to use approximate methods, which can be 
either analytical or numerical. 

The same problem holds for a quantum billiard, which 
is defined as a region of space where a point particle (for 
example an electron) is moving freely, not being allowed 
to escape. This is the quantum analog of a classical bil- 
liard where a particle only interacts when colliding elas- 
tically with the walls of the billiard, changing its trajec- 
tory. The Schrodinger equation in this case reduces to 
the Helmholtz equation used to describe the drum. 

The study of quantum billiards requires the calculation 
of a large number of energy levels, which are needed to 
extract the statistical distribution of the the level spacing 
Sn = En+i — En- On thcorctical grounds it is expected 
to obtain different statistical distributions if the corre- 
sponding classical billiard is either integrable or chaotic. 

Given the interest in calculating both the low and high 
parts of the spectrum, different approaches have been de- 
veloped over the years to obtain approximate solutions 
to the general problem. We may divide them into two 
categories: numerical and analytical methods. It is out- 
side of the scope of this paper to give a full account of 
all different approaches which have been used before for 
this problem; we rather prefer to mention a selected num- 
ber of works which we find representative and which can 
the source for additional reading. The reader interested 
in finding more details and references may find useful to 



check the bibliography of the papers that we mention; 
in particular, a very useful source of information on this 
problem is the beautiful paper by Kuttler and Sigillito 
which contains 142 references and which reports both 
numerical and analytical results. 

Among the numerical methods used to solve the 
Helmholtz equation we mention: 

• the method of particular solutions (also known as 
Point Matching Method) of Fox, Henrici and Moler 
01 (a recent improvement of this method is dis- 
cussed in ref. Q); 

• the method of fundamental solutions of Bogomolny 

i; 

• domain decomposition method of ref. Q (a recent 
application to isospectral drum is contained in [^); 

• the conformal mapping method (CMM) of Rob- 
nik 0; 

• the plane wave decomposition method (PWDM) of 
Heller 0]; 

• boundary integral methods (BIM) (see for example 

• the method of Vergini and Saraceno [ll| ; 

• the expansion method of Kaufman, Kostin and 
Schulten [H; 

• finite element methods (FEM) , see for example (isl . 
[l^ or the multigrid method of Heuveline [1^1 ; 

• collocation methods (in particular ref. contains 
a useful introduction to pseudospectral methods): 
as a special case in this category we mention the 
conformal collocation method [13, [3, which will 
be used in this paper; 

• the radial rescaling method of 
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We now come to the analytical approaches to the so- 
lution of the Helmholtz equation, which are essentially 
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based on perturbation theory. It appears in this case that 
the attempts to perform systematic and precise calcula- 
tions for the spectrum of arbitrary membranes are fewer 
and not as successfull as for the numerical methods. In 
particular Ref. [20j describes a boundary shape perturba- 
tion method, which is expressed in terms of the Green's 
function. Ref. [2]| describes the arbitrary perturbation 
of a circular domain which preserves the internal area of a 
membrane. This first order calculation is originally due 
to Rayleigh Otherperturbative approaches have 

been introduced in [H, Hi, [H, [H, [23|. Still based on 
perturbation theory is the A-cxpansion used by Molinari 
in [2^ to deal with regular polygonal membranes. None 
of the methods mentioned above has been applied to a 
systematic calculation of the spectrum, but rather of few 
selected levels. As a matter of fact, we are not aware of 
calculations in the literature for the whole spectrum of 
a drum, even slightly perturbed with respect to a solv- 
able problem. In a different context, the study of chaos 
in atomic nuclei, perturbation theory has also been used 
to study the statistical fluctuations of ground-state and 
binding energies [29| . 

Interestingly, it seems that the analytical methods have 
received considerably less attention than the numerical 
methods. The development of a systematic approach to 
perturbation theory for membranes (quantum billiards) 
of arbitrary shape is thus one the goals of the present 
paper. 

The paper is organized as follows: in section [II] we 
describe two numerical approaches to the solution of the 
Helmholz equation on arbitrary domains which are based 
on the use of conformal mapping (these are a new version 
of the CMM of ref. and the CCM of ref. [13); in sec- 
tion |TTT] we develop two different analytical approaches to 
the solution of the problem: in the first approach, which 
is nonperturbative, we prove two theorems which pro- 
vide, general expressions for the energies of the ground 
state or of an arbitrarily excited state; in the second ap- 
proach, which is perturbative, we obtain a systematic 
expansion for the energies of an arbitrary membrane. In 
section HVl wc apply the methods of sections HIl and HIH to 
several non trivial examples, which are helpful to illus- 
trate the power of the techniques developed in this paper. 
Finally in section |V] we draw our conclusions. 

All the numerical calculations have been carried out 
using Mathematica 7 [30| . 



II. NUMERICAL METHODS 



A. Conformal mapping method 

We briefly describe the conformal mapping method of 
Ref. The same notation of that paper is used here. 
V is the domain of the billiard whose boundary is given 
by an analytic curve in the (m, v) plane. The quantum 
states of a particle confined in V arc described by the 
hamiltonian 



H 



2m 



(1) 



with Dirichlct boundary conditions on V. In the follow- 
ing we will use fi^ /2m = 1. 

As explained in ref. a conformal tranformation can 
be used to map the original region V onto the unit disk: 



w = f{z) 



(2) 



where w = u + iv and z = x + iy. Actually, using Rie- 
mann mapping theorem, one may choose the mapping 
function f{z) so that V is mapped to an arbitrary sim- 
ply connected region of the plane, which we call fi. Here 
we choose to be a square of side 2. Upon performing 
the conformal mapping the operator in eq. ([1]) transforms 
as 



where we have defined 



dz 



(3) 



(4) 



We therefore obtain the eigenvalue equation 

A4,{x,y)+EJ:ix,yMx,y) = . (5) 

If we assume Dirichlet boundary conditions the eigen- 
fimctions iplx^y) may be conveniently expressed in the 
orthonormal basis of the eigenfunctions of the laplacian 
on the square Q as: 

1>n.,n„(a;,y)=sin(^(x + l)) sin(^(y+l)) ,(6) 

with UxjUy = 1, 2, . . . , 00. Von Neumann boundary con- 
ditions could be also studied by using the appropriate 
basis. We do not consider this case here. 

The eigenvalues corresponding to the wavefunctions of 
eq. ^ are 



In this section we describe two different numerical 
methods: in the first part we introduce and modify the 
conformal mapping method (CMM), which was origi- 
nally developed by Robnik [7|; in the second part we 
discuss the conformal collocation method (CCM), which 
has been developed by the author in a recent paper p7j . 



and we may write 



00 00 



(7) 
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By using the orthogonality of the <I> we obtain the equa- 
tions 

where 



dxdy $*/ „/ (a;, y)$„,,„„ (x, y)(10) 



In order to convert eq. ([9|) into a matrix equation we 
may identify a state in the box using a single integer 
index 



k = 



Nn 



y ' 



(11) 



where 1 < fc < N'^ and N"^ is the number of states to 
which we have limited the sums in eq. ([9|) ^. We may 
also invert this relation to obtain 



Tlx = 1 + 



k 



N + ri 



Uy = k — N 



N + 7] 



(12) 



where [a] is the integer part of a real number a and 77 — > 
0+. 

Using these conventions eq. ^ may now be cast in the 
form 



E 



^kk' 



1 

£k' 



^kk' 



Cfc = , 



(13) 



with 1 < fc < 7V2. 

The solution of this equation requires the calculation of 
a large number of integrals Y,kk' , which is the most time 
consuming part of the method, as explained in Ref. Q- 
Moreover, for large values of the indices in eq. (jlOp the in- 
tegrands are rapidly oscillating functions, which are more 
difficult to calculate numerically. 

The calculation however may be done quite efficiently 
if the conformal map f{z) is a polynomial or if it can be 
approximated with good precision by a polynomial, i.e. 
by the first few terms of its Taylor series around z = zq, 
with zq G D. In such a case we may express S(x, y) as 



00 00 

n— m— 



(14) 



and thus obtain 

00 00 

^n-j. ,ny — ^ ^ ^ ^ f^nniQur^n'^n Qnyji'^m (1^) 



n— m— 

where we have defined 



Qnmk^ J ax X sm(^— (x + l)J sm\^—(x + l)j 



^ Of course the indices in the equation span an infinite range, how- 
ever, the numerical solution of the equation requires to work with 
a finite number of it. 



The analytical calculation of these integrals may be 
done quite efhciently using the recurrence relations for 
the Qmnk given in appendix [A] Recurrence relations for 
the integrals which are relevant for the von Neumann 
boundary conditions, TZnmk, are also obtained. 



B. Conformal collocation method 

An alternative strategy for the solution of the 
Helmholtz equation is to use a collocation approach, as 
explained in Ref. [l3l- We briefly review this approach. 

As we have seen in the previous section the homoge- 
neous Helmholtz equation on a domain may be mapped 
to an inhogencous Helmholtz equation on a " simpler" do- 
main Q (which in this paper we assume to be a square 
of side 2). The resulting equation, eq. ([5]), may then be 
written in the equivalent form 



1 



T.{x,y) 



Ailj{x,y) = Eip{x,y) 



(16) 



The discretization of this equation proceeds as follows 
(more detail can be found in Ref. first, we intro- 

duce a set of functions, the little sine functions (LSF) of 



ref. [3l|, Sk{h,N,x), which are defined on x G (—1,1) 



and obey Dirichlet be at the borders (LSF correspond- 
ing to more general boundary conditions are obtained in 
ref. HI) : 



Skih,N,x) 



+ lix-kh)) 

2N 1 



cos 



1 sin(2^(x- fc/i)) 

((1 + w) U^ + kh)) \ ^^^^ 
coa{^{x + kh)) J 

For a given N (integer) there are iV — 1 of these functions, 
each peaked (with value 1) at a point Xk = kh {h = 2/N 
is the spacing of the uniform grid which these functions 
introduce) and vanishing at the remaining grid points 
Xj ~ jh, j ^ k. A function f{x) obeying Dirichlet be 
may be interpolated using the Sk{h, N, x) as 

/(^) 



N/2-1 

X f{xk)skih,N,x) 

k=-N/2+l 



(18) 



In a similar way we may derive twice this expression 
to obtain 



d'fix) 
dx"^ 



dx"^ 



d'^Skix) 



N/2-1 

k=-N/2+l 

N/2-1 N/2-1 

E E f(-^) 

k=-N/2+l j = -N/2+l 
N/2-1 N/2-1 

E E /(a^fc) ci'; s,(/i,iV,.T), (19) 

k=-N/2+l j=-N /2+1 



dx^ 



Sj{h, TV, x) 
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where in the last Hne we have introduced the matrix 

^(2) ^ cPsk(x) 



which provides a representation for the 



second derivative operator on the grid. Notice that these 
matrix elements arc known analytically. The discretiza- 
tion of eq. (|16p is now straightforward: for simplicity we 
consider the one-dimensional version of this equation and 
we write 



1 



I](a;) dx^ 



Skih,N,x) « -2^ 



cg^ s,ih,N,x()20) 



which provides an explicit representation of the operator 
O = on the uniform grid. 

If we want to generalize these results to two dimen- 
sions, we may consider the set of functions which are 
obtained from the direct product of the set of functions 
in the x and y directions, i.e. Sk{h, N, x)sk' {h, N, x). No- 
tice that we are assuming that the number of elements in 
the two directions is the same although different values 
could be considered if needed. Clearly the representation 
of an operator on the grid in terms of these functions 
will depend on four indices (the two integer indices cor- 
responding to the initial point on the grid and the two 
integers indices corresponding to the final point on the 
grid). However we may identify a specific point on the 
grid with a single integer 



N f N 

K^k' + — + {N-l){k+ — -l 



(21) 



where 1 < K < {N — 1)^. This relation is easily inverted 
to give: 



N 
N 



K 



N -1 + 7] 



k' = K~--{N-\) 



K 



N -1 + 7] 



(22) 
(23) 



where [a] is the integer part of a and 77 ^ 0"*" . 

Using these results one may easily generalize eqn. ([20)1 
to two dimensions as 



4? 4'/ + Skj c^',, s, {h, N, x)s,, (h, N, y) (24) 



from which we easily obtain the representation of O 
^ rA on the grid 



Okk' = - 



1 



^(2). , r A-2) 



(25) 



T.{xj,yj') 

We underline some useful features of this expression 



^ Eqns. II22I 1 and II23I I and the corresponding equations for the 
indices j and j' allow to express everything in terms of K and 
K', since k = k{K), k' = k'{K) and j = j{K') and / = j'{K'). 



,(2)> 



SkjC 



(2) 
k'j' 



• the matrix corresponding to 

the representation of the 2D laplacian operator on 
the uniform grid in a square of side L — 2; this ma- 
trix is "universal" , i.e. not specific to the particular 
problem one is solving, and sparse. 



• the matrix corresponding to 



is specific to 



^(xj,yji) 

the domain D since E is related to the conformal 
map; however this matrix is diagonal and therefore 
only — 1 elements need to be calculated; 

• the matrix O is a nonsymmetrical real square ma- 
trix whose eigenvalues are real. 

• no integrals need to be calculated, since E is simply 
evaluated on the grid points. In this case one docs 
not need to approximate the map with a polyno- 
mial in z. 

Because of the features above, the optimal computa- 
tional strategy consists of calculating the matrices for the 
discretized laplacian on fl for grids of different sizes. Typ- 
ically, this is the most time consuming task but it only 
needs to be done once, since these matrices are common 
to all problems and therefore they can be stored and 
used when needed without having to recalculate them. 
When a particular conformal map is chosen, one only 
needs to calculate the iV — 1 diagonal elements of E, 
which has a much more limited computational cost. The 
matrix representation for O is thus obtained quite ef- 
fectively and selected eigenvalues (eigenvectors) may be 
calculated rapidly. We will discuss some applications of 
this approach in Section HVl 



III. ANALYTICAL METHODS 

In this section we derive useful analytical approxima- 
tions for the energies of a particle confined in a two di- 
mensional region. In the first part we obtain a system- 
atic approximation to the ground state energy and wave 
function (or more in general to the lowest state in each 
symmetry class of the problem) : this expansion is really 
non-perturbativc and it is proved to converge to the ex- 
act values; in the second part we obtain the analytical 
expression for the spectrum of an arbitrary two dimen- 
sional region, worked out in perturbation theory up to 
third order. 



A. Nonperturbative formulas 

We will now derive a general formula for the ground 
state of a drum obeying Helmholtz equation ^ or the 
equivalent equation pB]) . 

As noticed before the conformal map has allowed to 
convert the original homogeneous Helmholtz equation, 
defined on a domain P, into an inhomogenous Helmholtz 
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equation, defined on a "simpler" domain fl, which in the 
present case we will assume to be a square of side 2. 

We niake few observations regarding the composite op- 
erator O = — -^(l y-^ A, which will help us to devise a suit- 
able approach to finding approximations to the states of 
eq.dlll): 

• the operator O is non-hcrmitian, 7^ O. This 
confirms the observation made in the previous sec- 
tion that the matrix representing the operator on 
the uniform grid is real and non-symmetrical. 

• the operator O may be cast in a manifestly hermi- 
tian form by considering its symmetrized form 



Qisym) ^ 1/^)1 



(26) 



In the framework of the Conformal Collocation 
Method (CCM) of the previous section, the matrix 
representing this operator on a uniform grid is sym- 
metrical and hcrmitian. Moreover its eigenvalues 
coincide with those of the unsymmetrized matrix. 
From now on we will always work with O^*^™) and 
therefore we will omit the superscript. It is inter- 
esting to notice that the form that we have found 
for O^^y"^) has been proposed by Zhu and Kroe- 
mer in [s^, discussing the problem of the connec- 
tion rules for effective mass wave functions across 
abrupt heterojunctions. A discussion of kinetic op- 
erators containing more general position dependent 
effective mass terms is discussed in (sj- 

• 51 is a two-dimensional domain on which the spec- 
trum (energies and wave functions) of the negative 
laplacian operator is known exactly, i.e. a square 
(rectangle) or a circle. 

• the eigenvalues of O are definite positive if we as- 
sume Dirichlet boundary conditions: the inverse 
operator O exists. For von Neumann boundary 
conditions a state with zero energy is present and 
one needs to define the invertible operator 



(27) 



where 77 —> 0+ is a positive infinitesimal parameter. 

• the spectrum of O is not bounded from above, i.e. 
there is an infinite number of states, of arbitrarily 
high energy; 

• the spectrum of the inverse operator is 
bounded from above, the largest eigenvalue being 
the reciprocal of the lowest eigenvalue of O (of 
course there is still an infinite number of states 
which become denser and denser as the eigenvalues 
approach zero). The spectrum of correspond- 
ing to von Neumann boundary conditions is also 
bounded from above, although its largest eigen- 
value tends to infinity as ry tends to zero. 



To simplify the notation in the following we will write 
Ojj as simply O, implicitly assuming that an infinitesimal 
positive shift rj is used if von Neumann be are used. 

We are now in position of stating the following Theo- 
rem: 

Theorem 1 Let 1^*0) be the exact ground state of the 
operator O and j^") be an arbitrary state with nonzero 
overlap with {"^o) ■ Then the lowest mode of the Helmholtz 
equation fil6]] is given by 



lim 

■n — ^00 



O' 



I*) • 



(28) 



C/(^) 



where Q-^ = \/S(-A + r/)"i\/S and E e 

is the conformal map which maps T) into More in 
general, if \^) is orthogonal to the first n states of eq. I116\) 
then the resulting state converges to the n + 1 excited 
state. 

The proof of this theorem is easy: for simplicity we 
limit to Dirichlet be and therefore set 77 = (the case of 
von Neumann be is obtained straightforwardly, working 
with a nonzero 77). 

Given an arbitrary state j^"), we may decompose it 
in terms of the eigenstates of O, \'^n), and write j^") = 

By applying the inverse operator to this state we obtain 



(29) 



where En are the eigenvalues of O. These eigenvalues 
are positive definite. En > 0, for all values of n. For 
this reason the component of j^I'o) in this new state has 
been amplified with respect to all the other components 
\^n>i)- Therefore repeated apphcations of the inverse 
operator will further suppress these components and the 
resulting state will converge to the true ground state of 
O. This concludes our demonstration. 

What we have described so far is actually an imple- 
mentation of the well known Power Method [s^ for an 
infinite dimensional matrix, with positive definite eigen- 
values which are bounded from above. 

There are several aspects which make the present im- 
plementation of the Power Method particularly appeal- 
ing: 

• the convergence rate of the algorithm depends on 
the ratio of Ei/ Eq, becoming slower when Eq and 
El are not well separated. An upper bound to this 
ratio is provided by the Payne-Polya- Weinberger 
conjecture [36l.[37t 



El 

Eo 



< 



El 
Eo 



disk 



jo,i 



2.539 



(30) 



where jo,i and ji^i are the first positive zeroes of 
the Bessel functions Jo{x) and Ji(x). This con- 
jecture has been proved in ref. [ssj . On the other 
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hand, it is easy to convince oneself that the ratio 
will get close to 1 for elongated domains 2?, where 
the longitudinal dimension is much larger then the 
transverse one. In these cases one may combine 
the method with a suitable dilatation in order to 
improve the convergence. 

• the ground state of the negative laplacian on Q pro- 
vides an initial guess of good quality (the typical 
implementation of the Power Method uses a ran- 
dom initial guess); 

• in the case of von Neumann be the initial guess 
should be chosen orthogonal to the zero energy 
state, which is trivial. 

• we may reasonably expect that the lowest modes 
of the negative laplacian on fl, which we call |n), 
dominate in |4'o); therefore we may leave the coeffi- 
cients of the expansion of the initial guess in terms 
of the |n) unspecified and use the variational prin- 
ciple to obtain them. In other words we may look 
for the coefficients which minimize the expectation 
value of O in the state generated with the Power 
Method. 

We now explore explicitly the application of the vari- 
ational principle to first order. For simplicity we assume 

to be a square. Let us choose an arbitrary state |x) 
which can be decomposed in terms of the eigenstates of 
the Laplacian on a box as 



ix> = E4"'ifc> 



(31) 



Notice that the state \k) = {k^, ky) is the direct product 
of the states corresponding to each direction. Similarly 
we call Ck the eigenvalue of the negative Laplacian on 
for this state, Cfc = ^ (fc^ -I- ky) . 

Working to first order we obtain an explicit approxi- 
mate expression for the ground state 



Using this result we may obtain the first order approx- 
imation to the ground state energy ^: 



(1) 



(35) 



where 



,(1) 



(1) 



(*J,"|0|*o 



oo 

E 

k,l.7n—0 



(0) (0) 



JO) m 



l.l' .7n,m'—0 



'-{m'\i:^/^\i'){i'\m) 



TO) 



(37) 



E'q^' given above provides an upper bound to the true 



energy Eq: 



En^^ > Eq. 



(38) 



The coefficients ch^ may therefore be chosen so that 
the expectation value of O is minimized. Therefore using 
a suitable number of coefficients in principle one may 
obtain arbitrarily accurate estimates for Eq. 

If we are willing to trade the accuracy for the simplic- 
ity, we may set c'q^ = 1 and c,^f>i = 0: in this case we 
obtain the approximation 



E. 



(1) 



EOO J_ 
1=0 ci 



|I]1/2|/)(^|I]1/2|0) 



En'=o^(o|sv2|0(;'|s|0(z|si/2|o) 



(39) 



which reduces to a quite simple formula when the internal 
sums are limited to the lowest state 



E, 



(1) 



(0|S|0) 



(40) 



In Section HVl we will consider few applications of the 
formulas given above. 

The generalization of these results to higher orders is 
trivial after observing that 



I* 



o- 



= ^ci")|fc) (41) 



14- 



(1) 



[6-1] [s^/^(-A)-ll]l/2] I;,) 



.(0) 



= y: —\k){k\^'^'\i)m' 



/2| 



6/ 

A; J,m— 

We may cast this result in the form 



where 



1^ 



^ JO) 



Lm,=0 ^' 



(32) 



(33) 



(34) 



^ The reader may notice that a variational bound on the ground 
state energy could also have been obtained calculating the ex- 
pectation value of O in |x): 



ixlx) 



Y^°° „(0) ■ 



and then minimizing this expression with respect to the c^^^ . 
This expression, however, has two unpleaseant features: first, it 
contains matrix elements of the operator S"-'/^, instead of E-*^/^; 
second, and most important, it contains the energies of internal 
states in the numerator, instead that in the denominator. Since 
in any practical application of these formulas a cutoff is used for 
the internal sums, this expression will be more sensitive to this 
cut and therefore less precise. 
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where 



oo (n—l) 



l.m=0 



(42) 



This equation provides a useful recurrence relation for 
the coefficients Cfc. 
For example we have 



oo (1) 

l,rn—0 
oo (0) 

/,m— 

where the completeness of states = 1 ^i^^ been 

used to obtain the final expression. 

The results that we have just obtained for the ground 
state may also be generalized for the excited states. As 
a matter of fact we may state the following Theorem: 

Theorem 2 Let I'i/if^) be an eigenstate of the operator 
O with degeneracy d and with an eigenvalue En- Let \^) 
he an arbitrary state with nonzero overlap with at least 
one of the \^\?) and let K be a real parameter for which 
\E„ 



A| <C 1. Then the state 



lim 

n — >oo 



O-A 



I*) 



(44) 



is a linear combination only of the d degenerate states 

The proof of this theorem is analogous to the proof that 
we have given earlier for the ground state. We just need 

to notice that the eigenvalues of the operator ( O — A 



are positive definite and bounded from above. Because 
of the condition \E„ — A| ^ 1 the repeated application 
of this operator to j^*) amplifies the components cor- 
responding to the eigenvalue En since l/(i?„ — A)^ is 

the largest eigenvalue of (^O — A^ . When the limit is 
taken only the components corresponding to this degen- 
erate eigenvalue survive, which proves the theorem. By 
choosing other d—1 linearly indipendents ansatz and re- 
peating the procedure, we obtain d linearly independent 
combinations of the degenerate states, which we can use 
to calculate the matrix elements of O. The diagonaliza- 
tion of this matrix provides the d eigenstates of interest. 
Notice that this strategy for calculating arbitrary eigen- 
values and eigenvectors of large matrices is discussed in 
refs. fHHSl. 



B. Perturbation theory 

While the previous approach based on the power 
method provides a systematic approximation to the 



ground state of a drum of arbitrary shape, an alterna- 
tive approach consists of applying perturbation theory 
to obtain the corrections to the energies and wave func- 
tions. 

Once again our starting point is the symmetrized op- 
erator 



(45) 



which reduces to the negative laplacian for S = 1, corre- 
sponding to a square of side L — 2. If we consider small 
deformations of this square and write 



S,, = 1 + 1](T 



(46) 



where cr = S — 1 and rj is a power counting parameter 
(which at the end of the calculation is set to 1), which is 
used to keep track of the different orders in powers of a. 
Our operator may now be expanded in powers of rj as 



O « Oo+vOi+if02+v^03 



(47) 



where the explicit form of the Oi may be worked out 
rather easily: 



Oo 



-A 



Oi = K-A) + (-A)a] 



(48) 
(49) 



62 = - [2a(-A)a + 3a2(-A) + 3(-A)a2] (50) 



63 = --[a\-A)a + a{-A)a^] 



(51) 



Within perturbation theory we may now calculate 
the corrections to the eigenvalues and eigenfunctions of 
our operator O, keeping in mind that it contains arbi- 
trary powers of rj. The standard scheme of Rayleigh- 
Schrodinger perturbation theory (RSPT) may be easily 
adapted to this problem and the contributions up to third 
order read: 



\{n\d,\k)\^ 



Cfc 



(52) 
(53) 

(54) 



|03|n)+2E 



{n\d2\k){k\di\n) 



+ EE 



n\Oi\m){m\di\k){k\di\n) 



{tn — ek){^n — Cm) 



|Oi|n)E 



{n\Oi\ky 
(e™ - Cfc)^ 



(55) 



Notice that the expressions corresponding to orders 77^ 
and rj^ receive contributions from distinct operators Oi. 



8 



We have calculated the perturbative corrections to the 
energy up to third order for the non-degenerate part of 
the spectrum: 



E^^^ = -en{n\(7\n 



(56) 

.^H^ (57) 



= -e^{n\a\nf + el{n\a\n)Y, 



{n\a\k) 



,\2 



UJ 1 



iel{n\(j\n) ^ 



k^n 



(n|g|fc)^ 



{n\a\k) {k\a\m) (m|cr|n) 



e"„ } } (58) 



Details of the derivation of these coefficients are given in 
the Appendix iBl 

The expressions above seem to suggest the presence 
of a term {—!)''' {n\(j\n)^ in the perturbative coefhcient 
of order k: assuming that this is true, one recovers a 
geometric series which produces a term en / {I + {n\a\n)) , 
which, for n = 0, is precisely the term obtained for the 
ground state energy working to lowest order, eq. (j40p . 

We may now discuss the extension of these results for 
degenerate states. Suppose that a state is d times degen- 
erate: in this case the perturbative expansion needs then 
to be carried out using the eigenstates of the dxd matrix 
ad of elements {ni\<7\nj)^ where \ni) (i = 1, . . . , d) is one 
of the degenerate states. With this modification the sums 
over internal states automatically exclude the degenerate 
states, given that the interaction does not mix them, and 
the perturbation expansion is well defined. Notice that if 
the matrix Wd is diagonal and proportional to the iden- 
tity matrix, the degeneracy of the levels will not be lifted 
by perturbation theory. We will discuss some examples 
of this in the Section IIVI considering the spectrum of a 
circular membrane. 

To illustrate how this works in practice we may con- 
sider the simplest case, in which the states are twice de- 
generate. Let us call \ni) and |n2) these states and ei and 
£2 their corresponding energies. If we work with pertur- 
bations of the square, the states of double degeneracy 
(which are also the most common states) are states with 
quantum numbers (nx,ny) and (jiy^nx). 

The eigenvalues of the matrix ad=2 are: 



£l,2 



\ [{ni\(T\ni) 



(?i2k|n2)±a2] (59) 



where 



62 = \j{{ni\G\ni) - (n2|CT|n2))' + 4(ni|a|n2)2 . (60) 

We clearly see that if {'ni\a\n2) = 0, then the £1^2 are 
degenerate (as anticipated). The eigenvectors of o'd=2 



Its eigenvectors are: 

\n1.2) = /3r''Vi)+M'''V2) , (61) 

where 




ni\a\ni) - (712 0- ^2 ) 



(62) 



^(2) ^ 1 /i + ^n,\a\n,) ~ ^n,\a\n2) ^^^^ 
V2 V 62 



^(2) ^ _l L_{nMn{)-{nMn2) ^^^^ 
V2 V 62 

We may now explicitly write down the different orders 
in PT. For example to first order we have 



£11 



[(7ii|a|7ii) + (712^1^2) ±62] . (66) 
To second order we have: 

^(2) _ , /- \2 , ,2 <>'^1,2| 



\a\k)^ 



2 

ve have: 

er^(^l,2|o■|?^l,2)^ + 4 X! 

k^n 

^ [(nl|cr|?^l) + (n2|cr|n2) ± 62]^ 

E ^ 1 M67) 

k^n 



We do not write explicitly the formula of the third 
order which is obtained in a similar fashion. 



C. Perturbation theory again 

Let us now try to obtain perturbative expressions for 
the energies using a different approach, which uses the 
theorems of eq. ([^5]) and (jH]). For simplicity we limit 
our considerations only to the ground state. 

Using eq. (pS)) we write 



Eq = lim 



I*) 



where 



0-1 = yrr^ (-^"') ^ 



(68) 



(69) 



and 1^') is an arbitrary state with non-zero overlap with 
the exact ground state. 

The perturbative series for Eq is of the form 



E, 



(PT) 



(70) 



fe=0 
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and therefore 



1 



(71) 



77=0 



to consider that the state j^") is arbitrary and therefore 
it may be chosen to be independent of 77. A suitable 
choice is to pick \'^) the ground state of the unperturbed 
hamiltonian, —A. 



In taking the derivatives of the expression (j68p we need 



Let us consider, for example, the first order term: 



E, 



(1) 



dE, 



(PT) 



drj 



lim 



m [o- 



ri=0 



m [o- 



(72) 



r;=0 



In order to evaluate this expression we only need to 
consider the terms containing the derivatives, since the 
remaining terms may already be evaluated at 77 = 0. For 
instance 



m[o-'T\^) 



E, 



(0) 



(73) 



Let us now consider O " us to order rj: 



o- 



« (-A)-" + 7? [-(-A)-" + (-A)-V(-A)-"+i 
+ ... + (-A)-"+V(-A)-i + (-A)-"- 



Therefore 



(n + l)(*|a|*) 



?)=0 



E, 



(0) 



(74) 



We thus obtain: 
E^^^ = lim [ni;o(*|(T|«') - (n+ l)So(^'|cr|*)] 

n — >-oo 

= -Eo{^\<j\^) , (75) 
which agrees with the formula obtained earlier. 

IV. APPLICATIONS 

In this section we consider several applications of the 
methods described in the previous sections and compare 
the results obtained with those in the literature. 

A. Unit circle 

The exact wave functions and energies of a circular bil- 
liard are known analytically and are expressed in terms of 
the Bessel functions and its zeroes. This example there- 
fore provides an ideal test both for our numerical meth- 
ods and for our analytical approximations (clearly O in 
this case is a square). Fig. [T] displays the unit circle 



which is obtained from a square of side 2 using a confor- 
mal map. The internal lines in the circle correspond to 
the mapping of a uniform grid on the square. Here the 
map has been approximated with the Taylor expansion 
to order z^^. The effect of this approximation is well 
under control, since the coefficients decay quickly: for 
instance, the coefficient of z^'^ is 1.6 x 10^^^. Wc may 
roughly expect that the smallest errors that we can ob- 
tain in our calculation stopping to this order are of this 
order of magnitude. Notice that using higher orders in 
the Taylor expansion makes sense only if one can han- 
dle the matrices of the size needed to achieve the extra 
precision made available by the Taylor expansion. 

We have applied both methods described in this paper 
to this problem, working to different orders. The results 
obtained for the first three levels are reported in Tables 
H] (CMM) andin (CCM). In the case of CCM is the 
number of elements of the basis of the square which have 
been used; in the case of the CCM iV^ is the number of in- 
ternal grid points used for the discretization. The reader 
may notice that the results obtained with the CMM for 
these states are more precise than those of the CCM for 
the same states: moreover the lowest error reached with 
the CMM is of the order of 10^". 

A second feature, which is common to both methods, 
is that the sequence of values obtained for a given eigen- 
value changing N is monotonously decreasing, as one 
would expect from the variational principle. This im- 
plies that the extrapolation of values corresponding to 
different N may be used to obtain more accurate results. 
Notice that this is not generally true for finite element 
methods, which do not describe exactly the border of the 
membrane. 

In Fig. [21 we have plotted the error ^ = ^^^^at — 1 for 
the first 3000 energy levels of the circle using either the 
CMM with N = 60 (thick black line) or the CCM with 

= 100 (thin red line). We need to make two observa- 
tions at this point: first, that the problem is symmetric 
with respect to both axes and therefore one could use 
this property to reduce the dimensions of the matrices to 
work with (for example, the even-even part of the spec- 
trum may be obtained in the CMM using a 900 x 900 ma- 
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trix); second, for a general problem the matrices obtained 
with the CMM are not sparse, while those obtained with 
the CCM are sparse. The sparsity of a matrix is a rather 
valuable property from a computational point of view 
since it allows to save computer memory. If we look at 
the two curves in the Figure we may confirm our previ- 
ous observation that the CMM provides better results for 
the lowest part of the spectrum: as a matter of fact the 
error for the ground state obtained with the CCM is of 
the order of lO"'^. This may be understood remember- 
ing that within the collocation method no integrals are 
calculated, since all expressions are simply evaluated on 
the grid. 

On the other hand the CCM results are obtained with 
a much larger matrix: for this reason we sec that the error 
grows more gently for the highly excited states compared 
to the CMM. 

We now consider the application of eq. psp : we may 
pick the coefficients c„ minimizing the expectation value 
in the equation, thus obtaining a precise upper bound 
to the energy of the ground state. In Fig. [3] we plot the 
error over the energy of the ground state calculated vari- 
ationally. A fixed cutoff Nint = 36 has been used for 
the internal sums; a varying cutoff N is used to limit the 
number of " final" states ■(/;„ (x) used in the calculation (re- 
member that the basis in two dimensions is obtained with 
the direct product of these states and therefore grows as 

In Table Hill we display the energies of the twenty lowest 
levels of the circular membrane obtained using perturba- 
tion theory to different orders (up to third order). In the 
first column we report the quantum numbers correspond- 
ing to the states of a square box of side 2; in the second 
column we report the degeneracy of the states; in the last 
column we display the exact results, obtained from the 
zeroes of the Bessel functions. The formulas obtained in 
this paper have been applied limiting the internal sum 



10"^ \ 



10-*- 



N 

FIG. 3: (color online) Error over the energy of the ground 
state of the circular membrane obtained using eq. p5|) and 
minimizing the coefficients. A'' is the cutoff in the basis (only 
the wave functions ^pi{x), . . ., ipMix) are used). The line is 

the fit Hf/t = 0.058 x lO"^''^'^. 



taking A'' = 20. 

There are several observations that we can make look- 
ing at this Table: first, that these results, although pcr- 
turbative, are rather precise, expecially for the lowest 
lying states; second, that we are recovering the correct 
degeneracies observed in the spectrum of the circle: in 
some cases, such as for the states (1,2) and (2,1) the 
the interaction term (a) is diagonal and therefore the de- 
generacy is not lifted by perturbation theory; in other 
cases, such as for the states (1,3) and (3,1) the inter- 
action is not diagonal and therefore perturbation theory 
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FIG. 4: (color online) Energy of the fundamental mode of the 
unit circle obtained using perturbation theory to order two 
(triangles) and three (diamonds) as a function of the number 
of states in the internal sums. The solid line is the exact 
result. The dotted line is the result of eq. (|40)l . 




FIG. 5: (color online) Deformation of the square of side 2 
produced by the map f{z) = z + /20. 



splits the levels However, the lowest of the splitted 
levels approaches the level corresponding to (2,2), thus 
approximately reproducing the degeneracy observed in 
the spectrum of the circle. Finally, for some states we ob- 
serve level crossing caused by the interaction: this is the 
case, for example, of the states (1,5) and (5, 1), which are 
splitted by the interaction, yielding a lower state which 
lyies below the degenerate levels (3,4) and (4,3). 

In Fig. m wc display the energy of the fundamental 
mode of the unit circle obtained using perturbation the- 
ory to order two (triangles) and three (diamonds) as a 
function of the number of states in the internal sums. 
The solid line is the exact result, while the dotted line is 
the result obtained using the simple formula of eq. (|40p . 



B. Deformation of the square 

We consider the conformal map f{z) = z + az"^, for 
I a I ^ 1, which produces a deformation of the square of 
side 2. FigO shows the membrane obtained for a = 1/20. 

In this case we have 



(76) 



We may easily calculate the matrix elements using the 
integrals earlier defined: 

{n\cj\m) = 4Q;Q„^m^i Suymy 

)(.77) 



* As we shall see later working with a general deformation, to 
leading order in PT the degeneration is lifted only for states (we 
limit to twice degenerate states) \nxny) where rix + ny is even. 



In particular: 

{n\a\n) 
for n ~ m and 



4a^ 



(78) 



32m,n, ((-l)-x+«. _ i) 
a ^ dn^ 



7r2 {ml - nl) 
64TOj,ny ((-1)'""+"=' +1) 



7r2 {ml ..y 



mf, — nf,y 



Q4m,n^ ((_!)"».+». ^ i) 
7r2 {ml - nlf 



(79) 



for n m. Let us discuss explicitly the case of states 
which are twice degenerate: |ni) = \nx,ny) and |n2) = 
\ny,nx) where ^ ny. 
For such states we have that the matrix elements read: 



(ni|cr|ni; 
{hi\a\n2) 



{n2\(T\n2) 
{n2\cr\ni) 



3 

, 



1 



9 2 



Since the matrix is diagonal and proportional to the 
2x2 identity matrix, the perturbation in this case does 
not separate the degenerate states and the first order 
correction in PT is the same as for a the non-degenerate 
states: 



-en-a^ ( 1 



9 9 



(80) 



We notice however that En'' goes like a^, so we need 
to calculate the contribution which comes from the 
second order correction. As we can see from the matrix 



12 



elements of a written above, to order a only nondiagonal 
terms contribute. Therefore, provided that ^ we 
have: 



{n\a\m) 



S2mxnx ((-1) 



1) 



7r2(m2 



o 



We may thus write the contribution to order in the 
second order PT: 



£;(2) 



2l4 

Fin,) 



where 



ml ((-1)™=^+"- 
2 [nl - 



1) 



(82) 



(83) 



Actually these formulas hold also for states of higher 
degeneracy, since also in this case the interaction matrix 
is diagonal and proportional to the identity matrix. 

The series appearing in this equation may be calcu- 
lated analytically (we do not report the expression here 
because of its length), instead we report the first three 
values are: 



F(l) = 
F(2) = 
F(3) = 



3072 
49152 

^4 



1024 

" 65536 ' 
5^2 



248832 746496 



-0.0164827 
i 0.0012288 
0.000325359 



We may also calculate the leading behaviour of F{nx) 
for Hx — > oo. In this case we find: 



lim F{nx) 



0.0317087 



3072 nt 



(84) 



Notice that in the limit of highly excited levels, Ux 
oo, we may use the asymptotic behavior of F[nx) written 
above to predict that the energy calculated to order 
behave as 



En 



1 + -a^ 



6 



- 1 



.(85) 



Provided that we consider states with n,„ ^ Ux we 



obtain the following limit 
En « e 



4 2 



(86) 



thus being decreased by a constant factor with respect to 
the energies of the box. 

In Table[IV]we display the quantity {En-E^°'')/a^ , for 
different values of a and for the first 50 states. The ener- 
gies in the deformed box have been obtained using CCM 
with TV = 60. In the last column we report the lead- 
ing perturbativc result, corresponding to the analytical 




FIG. 6: 



500 



1 



1000 1500 2000 2500 3000 
n 



^(CCM) 



for a = 1/100 for the first 3000 

levels. The numerical CCM results have been obtained with a 
grid corresponding to = 100. The inset plot is the blow-up 
of the first 200 levels. 



formula which has been obtained in this paper. Notice 
how well the numerical results approach the theoretical 
one even for moderate values of a. Only 15 of these cor- 
rections are positive, the remaining corresponding to a 
lowering of the frequencies. 

^(CCM) 



Fig. [5] displays the quantity ^ = 



1 - 



for 



a = 1/100 for the first 3000 levels. The numerical CCM 
results have been obtained with a grid corresponding to 
N = 100. The inset plot is the blow-up of the first 200 
levels. 

In Fig. [7] we display the first 50 levels of the deformed 
square as a function of 1/a. The solid line arc the precise 
numerical results obtained using CCM with = 60; the 
dashed lines are obtained using the analytical formulas 
that we have obtained. 



C. A general formula 

Let us therefore consider a general conformal map 

oo 

/(z) = z + r?^^z" , (87) 



where £>„ are real coefficients of similar strength and 77 is a 
small parameter. Notice that gi is related to a dilatation. 
In this case we have 



^ = Qn {{x + iy)"~^ + {x - ly)"-^) 

n=2 
00 

+ ^ gnQmix + iyT-^ix + iy)""-^ 



n.m— 2 



13 



200 



W=100 




1 / a 



FIG. 7: (color online) First 50 levels of the deformed square 
as a function of 1/q. The solid lines are the numerical results 
obtained with CCM with A'^ = 60; the dashed lines are the 
analytical perturbative results obtained in this paper. 



oo n— 1 , 



n=2 k=a 



(l + (-l)"-l-'=) 

OQ 



n,m— 2 



We focus on the term of order rj and write the corre- 
sponding matrix elements 



\a\m) w V^QnCn- 



n-i-m-j-num. 



+ 0\v'] , (89) 



n=2 



where 



E 

fe=0 



n — 1 



-1-fe 



(l + (-l)"-l-'=) 



(90) 



Notice that an analytical expression for the C is ob- 
tained quite easily using the results in the appendix [X] 
In particular the diagonal matrix elements read 



{n\a\n) 



VS 63 



(4 



48 nt 



Q5 



15 



09 



16{nx - ny){nx + riy) [ir^n'^ny - 90 {n'^ + n^)) 



'896 (90 I 



fi 2 , 4 4 



nlnl + nl) 



4 4 4 
71- «>y 



K - nlnl + n*y)) 32' 



45 



(91) 



We may now consider the off-diagonal matrix elements 
between two degenerate states. For the twice degenerate 
states \n) = |1,2) and |m) = |2, 1) we have 

(12|cr|21) = (21|cr|12) = . (92) 

In this case the degeneration of the levels is not lifted by 
the perturbation as already discussed earlier. The same 



occurrs for all the states where rix and Uy add to an odd 
integer. 

On the other hand, states where rix and ny add to an 
even integer do have non vanishing matrix elements. For 
example, for the states |n) = |1,3) and |m) = |3, 1) we 
have: 



27g5 , 63 871^^ 765 g9 297 3869775 - 42840^^ + 327r« gi3 
(13|a|31) = (31|a|13) = -— ^ + -g ^-^2 

135 (-1859 77516725 + 20659438807r4 - 16336327r8 + 256tt^'^) gn , , 

+ ^ + ■•■ (93) 



For these states the degeneracy is lifted by the pertur- 



bation. 
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where we have defined A 



as in Ref. 0. Since 



the maps f{z) and f{z) are related by a dilatation, one 
can calculate the eigenvalues corresponding to the second 
map and easily obtain the eigenvalues corresponding to 
the first map as: 



En — 



En 



cos^p 



(96) 



We will therefore use the second map and then use the 
relation above to obtain the eigenvalues corresponding to 
the first map. We thus have: 



1 + 4Ar cos ( 



4A^r 



2„2 



(97) 



and 



FIG. 8: (color online) Robnik's billiard for A = 1/2. 



a = AXr cos 9 + 4A 



(98) 



In Table |V| we have written the analytical expressions 
for the first order perturbative corrections (divided by the 
energy of the square of corresponding quantum numbers 
and changed of sign) generated by the conformal map 
written above, for the first 20 energy levels. The reader 
may notice that only odd powers in the map contribute to 
this order. Moreover, the degeneration is lifted only for 
the states {|1, 3), |3, 1)}, {|2, 4), |4, 2)} and {|1, 5), |5, 1)} 
(for which as mentioned before the sum of the two quan- 
tum numbers is even), although the splittings depend 
only on some of the coefficients (£'5,^9, . . .). 



D. Robnik's billiards 

Out next application is to Robnik's billiard, which are 
obtained using the conformal map of the circle Q 



Az + Bz'^ 



(94) 



where A = cosp and B ~ ^ sinp. The parametrization 
of A and B in terms of p describes a family of billiards 
of constant area A = ir. We write the map as 



f{z) ~ cosp [z + Az^l = cosp /(z) 



(95) 



Notice that having chosen to work with the second map 
rather than with the first one, eliminates a constant term 
a and therefore simplifies the calculation. 

In this case we are working with the orthonormal basis 
of the circle: 



4'kns{r,9) = RknJkilknr) X 



cos k9 
sin k6 



s = 1 
s = 2 



(99) 



The normalization constants are: 



Ron — 



1 



/tt Jo(7o,i) 
^/2 



, fc = 
, fc > 0, 



(100) 
(101) 



where Jk is the Bessel function of order k and "fkn its 
zeroes. We have followed almost entirely the notation of 
Ref. 0, apart of introducing an index s for the degen- 
eracy: all states with A: > are twice degenerate, states 
with fc arc simple. 

The application of our perturbative formula requires 
the calculation of the matrix elements of a. Thus we 
need to calculate: 



/•I /•27r 

{kns\r^\k'n's') ^ / dr r''cl)kns{r,0)cbk'n's'{r,0) 

Jo JQ 

= TT[Sk.k' + hiy+^Sk.-k']Sss'RknsRkn's f drr^ Jk{jknr) Jkilkn' r) (102) 

Jo 

{kns\r COS 9\k'n's') = dr dcj) cos 0(j3knsir,6)(j3k'n's'ir, 6) 
Jo Jo 
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3kk' s^ss' Rkns 



Rk' 



drr^Jkijknr)Jkhkn'r) 



(103) 



where 

gkk's = [{5k-k\l + 5k-k' -l) 

+ (-i)''+Mw,i + 4+fe',-i)] . (104) 

Using the integrals above we may obtain the matrix 
elements of cr: 

{kns\a\k' n s) = \{kns\a''^^\k'ns)+\'^{kns\a^^'^\k'n's) 
= 4X{kns\r cos 6\k'n' s) 
+ 4:\^{kns\r^\kn's) . (105) 

Notice that the diagonal matrix elements are of or- 
der A^, whereas the off-diagonal matrix elements start 
at order A. This means that our pcrturbative (in the 
power counting parameter rj) formulas contain contribu- 
tions with different orders in A. For example the first 
order formula is of order A^, just as the second order per- 
turbative formula, which contains squares of off-diagonal 
matrix elements. On the other hand, the third order 
perturbative formula is of order A^, since the last sum, 
which involves the product of three non-diagonal matrix 
elements vanishes when the external states are taken to 
be equal. For this reason, the general results that we 
have obtained in this paper allow us to make a rigorous 
perturbative analysis of this problem only to order A^. 

We may therefore write 

Ekus=E^^i+X' Eil+0[X^] (106) 

where 

Ekns - Ikn 

^/ ^/ 'kn 'k'n' 

and "Y^'f.! „/ means the sum obtained excluding the values 
{k',n') = (fc,n). 

Even before doing any numerical calculation, we may 
look at the form of the matrix elements of a to under- 
stand how the degeneracy of the level of the circle is 
lifted. Wc observe that the s dependent component of 
gk,k'.s is nonzero only when (fc, k') is (0, 1) or (1, 0). This 
means that only the levels with principal quantum num- 
ber k ~ 1 have the degeneracy lifted. A splitting of the 
levels corresponding to fc > 1 may only occurr through 
higher order perturbative contributions (to order A* and 
higher, and it is therefore highly suppressed. 

In Table IVII we display the first 40 eigenvalues of a 
Robnik billiard with A = 1/100 and A = 1/20, calculated 
with the perturbative formula (|106p to order A° and A^ 
and with the CCM with N = 100. We may observe that 



the splitting of the levels corresponding to fc = 1 is repro- 
duced quite well by our formula; in particular the results 
obtained for A = 1/100 with second order perturbation 
theory reproduce amazingly well the numerical CCM re- 
sults. For A = 1/20, wc observe a small splitting of some 
of the levels corresponding to fc > 1, such as for example 
for the states with quantum numbers (2,1). As men- 
tioned before, these splittings come from higher order 
perturbative corrections, which we are not considering 
here. 

In Fig. [5] we display the splitting of the first two ex- 
cited levels of the Robnik's billiard, corresponding to the 
quantum numbers (fc,n) = (1,1), as a function of the 
parameter A. The solid line is the perturbative result to 
order A^, calculated in this paper, which corresponds to 
^£;(PT) ^ 29.3639418A2 + 58.7278836A4 (observe that 
this expression contains a power A* because of extra 
power of A coming from the rescaling factor 1 / cos^ p in 
the energies); the pluses are the precise numerical re- 
sults obtained using CCM with N = 100. Finally, for 
A 0.5, the splitting is seen to behave linearly in A: the 
dashed line represent the linear fit of the last four points, 
AE'^f'^') = 0.5095 4- 4.90231A. A resummation of the per- 
turbative series (assuming that it is possible to calculate 
it to arbitrary orders), should therefore be able to repro- 
duce this leading asymptotic behaviour. It would also 
be interesting to study the radius of convergence of the 
perturbative series and see if it has any relation with the 
insurgence of quantum chaos. 



E. Regular polygonal billiard 

Another interesting application of the techniques dis- 
cussed in this paper is to the calculation of the frequen- 
cies of a regular polygonal billiard. This problem has 
been attacked in a series of papers, using both numerical 
and analytical approaches. 

In ref. [13 Liboff has studied the polygonal quantum 
billiard problem, proving that the first excited state of a 
regular polygon is doubly degenerate. In ref. [i^ Liboff 
and Greenberg have studied the hexagon quantum bil- 
liard, constructing a subset of eigenfunctions and eigen- 
values, while in ref. [i^l Liboff has considered classical 
and quantum wedge billiards. 

Molinari has obtained analytical expressions for the 
frequency of the ground and first excited states using a 
perturbative approach, which he calls "lambda expan- 
sion" [2^. Cureton and Kuttler [i^ have used a numer- 
ical approach (which is essentially the Conformal Map- 
ping Method of 0), to obtain quite precise values for 
selected eigenvalues of regular polygons and of the fig- 
ures obtained by the dissection. More recently, Grinfield 
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where /o = 1 and 
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(109) 



for fc > 1. 

As explained in the previous section, it is convenient 
to work with a rescaled map defined through the relation 
f{z) = Cj\ff{z) and then easily relate the eigenvalues 
corresponding to the map /(z) to those corresponding to 
the map f{z) by means of the equation 



Er, 



En 



(110) 



We thus obtain 



FIG. 9: (color online) Splitting of the levels corresponding to 
(k,n) = (1,1). The solid line is the perturbative result to 
order A'^; the pluses are the numerical results obtained with 
COM with iV = 100. The dashed hue is the linear fit of the 
four last numerical values AE'^-'^'*' = 0.5095 + 4.90231A. 



and Strang [4y| have obtained an analytical formula for 
the simple states of regular polygonal with M sides as a 
function of \/M. 

Other works have considered specific polygons, and ap- 
plied different methods to solve the Hclmholtz equation 
on these domains: for example, Betckc and Trefcthen 
have applied the method of particular solutions in rcf. 
to obtain very precise values for the first three eigenval- 
ues of a regular decagon; Lijnen, Chibotaru and Ceule- 
mans, used a radial rescaling approach in rcf. to ob- 
tain the ten lowest eigenvalues of regular polygons with 
Af ~ 3,..., 6 sides; finally, Guidotti and Lambers [i^l 
have considered a variant of the method of particular so- 
lution, and applied to a regular polygon with 128 sides 
(comparing their findings with those predicted by the 
formula of Grinfield and Strang). 

We will now treat this problem with the systematic 
perturbative approach developed in this paper and then 
perform a comparison with the works mentioned above. 

Our starting point is the conformal map 



!{z) = Cm 



ds 



(1 _ s^f)2/^f ' 



(107) 



which maps the unit circle into a regular polygon of N 
sides inscribed in the circle [1^, [1^ . Following Molinari 

we have defined C_\f = r ( i + 1 / aT) r^(i'~2 /j\f) have 
expanded the map as 



fc=0 



E = ^^(AAfc + l)(AAj + l)A./,r^(^+J)e*-^(^-^)^ 



fc=0 j=0 



and 



oo oo 



a = -l + ^^(AAfc + l)(AAj + l)/fe/, 



k=0 j=0 



Let us calculate explicitly the matrix elements 
{kns\a\k'n' s) . We define 



7-(i) 

k,k' ,/,m 



COS k9 cosk'9 cos JV{1 — m)9 d9 



[^k-k'+Afil-Tn.) + 5k+k'+N{l-m) 



j(2) 
k^k' ^Lni 



Sk~k'~Af{l~m) + Sk+k'-Af{l~m)] (HI) 

/ sinfc^ sink' 9 cos U [l - m)6 dO 
Jo 

= 2 [^k-k'+Af(l-m) — Sk+k'+Af{l-m) 

+ 5k-k'-J\f{l-m) — 5k+k'-J\f{l-m)\ (112) 

which can be written as 



-^k,k',l,m — -^k,k',Lm ' I ^) -^k,k',Lm 

= [^k-k'+Af{l-m) + ^k-k'-Af{l-7n)] 

+ (^l)''^"^'^ [Sk+k'+^r{l-m) + Sk+k'-ATH-m)] 



We also define 

\Jkn^k'n' ,M ,l^m - 



r^^'+"'^+^Jk{lknr)Jk'{lk'n'r)dr 



(108) and thus we may write the matrix element of a as: 
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{kns\a\k' n' s') = ~5kk'5nn'5ss' + (-^^ + l)(A/"m + l)flfmRknRk'n'Sss'llfl>,i,^Jkn,k'n',Ar,Lm ■ (113) 

/=0 "1=0 

I 

We focus here on the calculation of the first order cor- We may now explicitly write the diagonal matrix ele- 
rection to the eigenvalues: ment of a: 

eH = ~-/L{kns\a\kns) . (114) 
I 



OO 

{kns\a\kns) = + 1)^ f?RLli%^iJkn,knM,i,i 

1=1 

oo oo 

+ (-1)'+' E E (-^^ + + ^)flfrnRL^i!'lunJkn,knM,l,m , (115) 



where the last term, when non-zero, breaks the degen- 
eracy of the states with s = 1 and s = 2. If we look 
at the form of ij^'l. ^ ^ we see that it does not vanish 
only when the condition 2k = zLAf{l — m) is met. We 
need to distinguish between regular polygons with even 
or odd numbers of sides: for Af even, we may fulfill this 
condition for k = Af/2 and / = m ± 1, for fc = A/" and 



I = m zt 2, and so on; for odd, this condition is met 
for k ~ M and Z = m ± 2, for fc = 2I\f and I = nizLA, and 
so on. Thus, perturbation theory to first order predicts 
a peculiar pattern of lifting of the degeneracy of levels. 

Notice that the radial integrals can also be calculated 
explicitly: 



T{k + l)^{2k + Af{l + m) + 2) 
X 2F3[k + -,k+ — + — + l;k + l,2k + l,k+ — + — + 2-^-llA . 



We can therefore write the energy calculated up to first 
order in perturbation theory as 



^kns 



'-■AT 



[1 — {kns\a\kns)] 



(116) 



Let us focus for a moment on the perturbative result 

2 

obtained at order zero, Ekns ~ t^- We may expand the 
factor 1/C^ around A/" = oo, obtaining the energies 



pi 



1 



27r2 4C(3) 147r'' 



3AA2 AA3 45Ar4 



which reproduce the formula obtained by Grinfield and 
Strang in [iHl up to order l/A/"^ (notice however that 
their formula is limited to the simple states). The coef- 
ficient of order XjM'^ is different, although very similar 



numerically: we obtain 4^(3) ~ 4.81, while ref. [4Q| quotes 
7r'^/6 « 5.2 ^ . On the other hand, if we look now at the 
first order correction, we clearly see that it is not analyt- 
ical at A/" = c» and therefore that it cannot be expanded 
in powers of 1 jM . 

In Fig.[Tni[II]and[Tl]we have plotted the first 40 levels 
of the octagon, nonagon and decagon billiards, obtained 
using zero and first order perturbation theory (circles and 
triangle respectively) and using CCM with a grid corre- 
sponding to iV = 100 (pluses). The boxes in the plots 
enclose the levels where first order PT predicts that de- 



The value found by Grinfield and Strang is based on a conjecture 
based on numerical evidence. 
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FIG. 10: (color online) First 40 levels of a octagon billiard. 
The circles (triangles) are the results obtained with zero (first) 
order perturbation theory; the pluses are the exact (numeri- 
cal) results obtained using CCM with A'^ = 100. The boxes 
enclose the levels where first order perturbation theory pre- 
dicts that the degeneracy of levels is lifted. The corresponding 
quantum numbers are displayed. 



generacy is lifteid: as we have already (discussed before, 
this occurrs at k = Af /2 (and multiple integers) for regu- 
lar polygons with even number of sides and at fc = A/" 
(and multiple integers) for regular polygons with odd 
number of sides. The numerical values corresponding 
to these plots are reported in Tabic IVlTl we like to un- 
derline that in some cases the CCM may provide tiny 
splittings for states which are degenerate. This occurrs 
because the discretization may explicitly break the sym- 
metry of the billiard, although this splitting is bound to 
decrease as the grid gets finer and finer. Taking into 
account this fact, and looking at the Table we may spec- 
ulate that the states lying below the {Af/2, 1) states for 
even billiards and (A/", 1) for odd billiards are degenerate 
(remember that Liboff's theorem of ref. [i^l tells us that 
the first excited state of the regular polygonal billiard is 
degenerate). A rigorous proof of this conjecture would 
require the calculation of higher order perturbative con- 
tributions, which we plan to consider in the future. 



FIG. 11: (color online) Same as for Fig. 1101 for the nonagon 
billiard. 
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FIG. 12: (color online) Same as for Fig. [TO] for the decagon 
billiard. 



F. Hamiltonians with position dependent effective 

mass 



As we have mentioned in section IIII Al the sym- 
metrized form of the density dependent negative Lapla- 
cian operator that we have proposed in this paper had 
already been introduced in a paper by Zhu and Kroemer, 
ref. [33| , discussing a different problem. 

In recent years there has been lot of interest in dis- 
cussing the solutions of quantum mechanical models con- 
taining density dependent effective masses, which are 
common in the description of condensed matter systems. 
However most of the efforts have been put into con- 
structing exactly and quasi-exactly solvable potentials for 
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the Schrodinger equation with position dependent mass 
terms [i^. To the best of our knowledge perturbation 
theory has not been used as a tool to construct the solu- 
tions to this problem; since the results obtained earlier in 
our paper apply straightforwardly to the present case we 
now extend the approach of section llll BI to Hamiltonians 
of the form: 



H 



1 



2m 



1 

7^ 



V 



(117) 



Notice that the stationary Schrodinger equation asso- 
ciated to this operator may be solved numerical using 
the CCM of section III Bl with the difference that now 
E is a physical density, not related to a conformal map- 
ping. Moreover, the matrix representing the operator 
H on the uniform grid, now contains also a diagonal 
term corresponding to the discretization of the poten- 
tial V. For problems defined on whole space, a suitable 
strategy corresponds to solve the SSE in a d-dimensional 
cube {d being the dimension of the problem), whose di- 
mensions have been determined variationally, as done in 
refs. [m, . The only difference in this case is the pres- 
ence of the spatially dependent density S. 

An alternative way of solving the problem consists of 
resorting to perturbation theory. We may cast this hamil- 
tonian into the form 



H 



1 



2m 



A + Vb 



1 



V 



, (118) 



where Vb is a potential for which the "unperturbed" 
Hamiltonian Hq = ^— -I- Vq] is exactly solvable. 



The analysis carried out in section llll BI mav now be eas- 
ily generalized, by expressing Y, = 1 + rja as before and 

letting (V-^) ^ 1] (v - -r^] . Thus we have 



l+))cr 



where 
Hi = 

H2 = 

= 



+ V-Vo 



H. 



(a) 



H. 



(119) 

(120) 
(121) 



3 

^16 
Vo<J^ 



where the operators are those which contain Hq and 
are equivalent to the one obtained before. Working up 
to second order and calling e„ and \n) the eigcnenergies 
and eigenstates of we find 
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a' He, + Hna- 



(122) 



'en{n\a\n) + {n\V -Vo\n) 



(123) 



^ e,,{n\a\nr + elY.- 



n\a\ky 



n|i7i°-'|fc)(fc|ff}°V> 



E 



(n|gf^|fc)^ 



(124) 



k^n 

Let us consider the term: 



k^n 



-Y^^lJ-±l^n\a\k){k\V ~V,\n) 

, , ^k 
k^n 

- {{n\a{V - Vo)\n) - {n\a\n){n\V - Vo\n)) 



2E 



Therefore we have 



{n\a\k){k\V-Vo\n) 



(125) 



e„(n|cr|n)^ -|- el 



k^n 



{n\a\k)^ 



{{n\a{V - Vo)\n) ~ (n|a|n)(n|y - V^\n)) 

'^y^^^{n\o\k){k\V -Y^\n) 
f-f £« - efe 

k^n 

(n\V-V^\kY 



E' 



efe 



(126) 



Likewise we may find the perturbative corrections to 
the eigenstates of the Hamiltonian; for example, we may 
consider the special case where V = Vb/\/S, and there- 
fore obtain to first order in perturbation theory: 



fife 



2 ^r. 
k^n 



Cfe 



\k){k\a\ 



(127) 



We may apply this expression to the calculation of the 
uncertainties over position and momentum: 



Ax„ « Aa;l°) 



_ ?7 \p e» + £fc {k\(j\n) 

(a.!"))' 

X {l^n\x^\k) - {n\x\n)(n\x\k))) (128) 



Ap„ « Ap(") l-^E?^ 



e/c {k\a\n) 



X {{n\p^\k) ~ {n\p\n){n\p\k))) . (129) 
We thus obtain the product of the uncertainties 
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2 ^ e„ - Cfe 



[{n\x^\k) — {n\x\n) {n\x\k)) [{n\p'^\k) — {n\p\n) {n\p\k)) 



Ax 



(0) 



Ap! 



(0) 



(130) 



which for the simple harmonie oscillator simplifies to 

Ax„Ap„ « Axi°)Ap(,°) (l + 0[772]) . (131) 



V. CONCLUSIONS 

Our paper contains strong numerical and analytical re- 
sults for the solution of the Helmholtz equation on gen- 
eral domains. We briefly state the novelties and advan- 
tages of the approaches outlined in this work: 

• The extension of the CMM of Robnik to use 
the orthonormal basis of a square allows to express 
the integrals (for polynomial conformal mapping or 
mapping which can be well approximated by poly- 
nomials) analytically. This clearly reduces the com- 
putation times and improves the precision. 

• The computational advantages of the COM of 



Amore 171 are made clear: the collocation ma- 



trix is made of two parts, a "universal" matrix, 
corresponding to the discretization of the negative 
Laplacian, which can be calculated once and for 
all for a given grid, and a specific "shape depen- 
dent" matrix, which is diagonal and therefore may 
be calculated quite efhcicntly. In this approach no 
integrals are ever needed to be calculated. 

• We prove two theorems which provide explicit for- 
mulas for the ground and excited states of a quan- 
tum billiard. These theorems essentially implement 
the very well known "power method" for finite di- 
mensional matrices, to operators in an infinite di- 
mensional Hilbert space. The variational theorem 
also allows to obtain upper bounds to the ground 
state energy. 

• We have formulated a systematic perturbation the- 
ory for membranes obtained slightly deforming a 
square or circular membrane. Our perturbation 
scheme is essentially the Rayleigh-Schrodingcr per- 
turbation theory with a specific form of the in- 
teraction operator. The perturbative results dis- 
cussed in this paper are in most cases very precise 
and provide a clear and simple explanation of the 
mechanism responsible for the lifting of degenera- 
cies. For the case of the small deformations of the 



square, we have obtained an analytical formula for 
the whole spectrum. A possible application of our 
shape perturbation theory could be to the calcu- 
lation of Casimir energies due to deformations of 
membranes. In a recent interesting paper, ref. [50t . 
for example, the point matching method has been 
used to numerically calculate Casimir energies for 
perfect-conductor waveguide of arbitrary section. 
Our perturbative approach could be useful for an- 
alytical calculations. 

• The shape perturbation theory of the present pa- 
per could be applied to study the statistical dis- 
tribution of the energy level spacings for quantum 
billiards obtained from small perturbations of the 
square or of the circle. In this case the study 
would provide an alternative tool to purely numer- 
ical methods or to asymptotic methods. 

• We have provided a way to extract (at least in prin- 
ciple) the perturbative coefhcients of the energy of a 
given level to a specific order. It would be certainly 
interesting to study the convergence properties of 
the perturbative series corresponding to quantum 
billiards of different shape and see if a relation be- 
tween the chaotic behaviour of a billiard and the 
perturbative series could be established. 

• The analytical results contained in the paper are 
not limited to two dimensional membranes, but can 
be used for problems of higher dimensions; more- 
over they can also be applied to describe the spec- 
trum of arbitrary inhomogeneous membranes. A 
numerical study of this problem using the CCM 
has been already performed in [l^ ; we plan to look 
at this problem with the analytical techniques of 
the present paper in the near future. Here we 
have applied our results to Hamiltonians containing 
position dependent effective masses, extending the 
perturbative approach developed for the Helmholtz 
equation. 

We believe that several developments may stem out of 
the results contained in this paper: among these, it would 
certainly be worth looking at possible non-perturbative 
extensions of our shape perturbation method, which we 
plan to look at in the near future. Similarly, we plan to 
carry out an extension of the perturbative calculation of 
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this paper to higher orders. A further interesting apphca- 
tion could be to the inverse problem: in a classical paper, 
Kac [m posed the question whether one can hear the 
sound of a drum, meaning by this if the complete knowl- 
edge of the spectrum of a drum is sufficient to determine 
its shape. We now know that the answer to this question 
is negative, since in 1992 Gordon, Webb and Wolpert 
[5^ found an example of two two-dimensional not iso- 
metric drums which are isospectral. On the other hand, 
the perturbative formulas obtained in this paper might 
be used to find the shape of a drum which is obtained 
from a small deformation of a square or of a circle, know- 
ing a number of frequencies. Finally the last application 
that we wish to mention is to the vibration of drums 
with fractal boundaries: refs.fs^, [13, [EE], for example, 
consider the vibration of a membrane whose border is 
the Koch snowflake. In particular Banjai [5^ uses a con- 
formal transformation to map the snowflake into the unit 
circle and then performs a spectral collocation to obtain 
numerical results. Our perturbative approach, given the 
conformal map, is directly applicable to this problem and 
could maybe provide further insight in this problem. 
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APPENDIX A: RECURRENCE RELATIONS 

We introduce the definitions: 

+ (Ai) 

Xn{x) EE ^+0(„-l)cOs(^(x+l)) , (A2) 

where 6{n — 1) = 1 for n > 1. 

Notice that these functions obey Dirichlet and von 
Neumann boundary conditions respectively and that 
they are orthonormal: 

-1-1 r+l 

(t)n{.x)4>m{x)dx = I Xn{x)Xm{x)dx = 5nm i.^'i) 
-1 J-1 

We consider the integrals: 

Qnrnk = dx 4>n{x)4>m{x) (A4) 

r + l 

Unmk = y dx X^ Xn{x)Xm{x) . (A5) 

With simple algebra we are able to obtain recurrence 
relations. First we consider the case m ^ n: 



mn TT 



w? -\- ri^ TT^ 
"(fc-f-l)(fc-f 2)T 



Qnmk+2 



(fc + l)(fc + 2) 2 



■Ti-nmk+2 



(A6) 



,2 I „2 



_ (-l)" + (-l)"'+" m- + n- 



k + l (fc + l)(fc + 2) 4 

^2 



mn TT 



(fc + l)(fc + 2) 2 



Q 



(A7) 



These equations may be solved to yield the recurrence 
relations: 



Q 



nrnk-\-2 



8(fc + 2)mn ((-1)'= + (-1)™+") 

(^jj^2 _ 7^2^2 

4(fc + l)(/c + 2) {m^ + n^) 

r- ; ^"2 Qnrnk 

8{k + l){k + 2)mn 

— - —2 Ti-nmk 

4(fc + 2) (m2 + n2) ((-1)'= + (-1)'"+") 

^2 ^jy^2 _ j7,2^2 

8{k + l){k + 2)mn 

— ; —2 Qnrnk 

4{k + l){k + 2) (m2 + n2) 



(A8) 
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For m = n we have 



Qmnk ~t~ "^nmk 



fc + 1 



(A9) 



(AlO) 



and 



Qnrnk 



k , 1 A ^2 



2(fc3 + 6fc2 + ll/c + 6) 

,r2r,2 

Qnrnk+2 (All) 



F + 3fc + 2 

which is easily solved to give 

(-1)'^ + 1 k^ + 2,k + 2 



Q 



nrtik-\-2 ' 



2(A: + 3) 7r2n2 



Q„,„fc (A12) 



After calculating explicitly Q„„iO, Qnmi, Tinmo and 
T^nmi one may then obtain the coefficients of higher order 
by simply applying the recurrence relations. 

For instance: 

QnmO ~ ^nm 
Qnnl = , m = 71 

8mn((-l)'"+" - 1) 



Qnnil 
Qnn2 
Qnm2 



^2 ^j^2 „ j^2^ 

1 2 

, m = n 

16mri((-l)™+" + 1) 



m ^ n 



7r2 (rrfi — 'n?) 
Qnnz = , m = n 
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APPENDIX B: PERTURBATION THEORY 

In this section we give some detail on the calculation of 
the perturbativc coefficients for the energy. In the follow- 
ing we introduce the definition a;„fc = e„ — e^. We may 
use the expressions given in eans. ((i5|) . ([i51) . ([5n)l and ([5T|) 
inside the general formulas of cans. ([5^ . ([55]) . ([51)1 and ([55]) 
to obtain the analytical expressions which are reported 
in the paper. 

We start with the second order term, since the zero 
and first order terms are straightforward. In this case we 
obtain: 



£;(2) 



^(n|(7|0e;+3e„(n|f72|n) 



{n\a\kY 



4 f-^ LUnk 



(Bl) 



Notice that the first sum in this expression is unre- 
stricted, since I can take all possible value, including 
I = n. It is convenient to express the unrestricted sums 
isolating the term corresponding to / = n, i.e. writing 



Keeping in mind this observation we may rewrite the 
first term of the equation above as 
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-en{n\a^\n) + -e^(n\a\n)^ + - ^(:n\a\kftu 



Similarly we may simplify the second term: 
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Combining the two terms we obtain the final result: 
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The third order contribution reads: 
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We may simplify this expression by arranging term with equal number of sums together. Keeping in mind this 
simple observation we consider the terms in the above expression which contain two sums: 
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-I- Cfe) ek{€n + e,n) , (e„ -I- efc)(en + em)(em + Efc) 



^nk^nm 



4E E("Hfc)(fc|'^HH^i' 



k^n my^n 

We start with the first term in the square parenthesis 



2-3e„ ( — - 


. ' ) 


+ 4 ^" 






'-^nk^nm _ 



2surns 



-y {{n\a^\n) - 2{n\(j^\n){n\(j\n) + {n\a\nf) 



(B4) 



(B5) 
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We now come to the second and term terms in the square parethesis: a simple relabeling of the summed indices 
shows us that these terms are equal and therefore: 



II 



2sums 2 



EE' 



{n\a\k) {k\a\m) {m\(7\n) 



^ '^^-^[{k\a'\n)~{k\a\n){n\a\n)] 



2 f-r' (^nk 
k^n 



(B6) 



Finally we write the last term as: 



/// 



-?.EE 

k^7i m^n 



{n\a\k) {k\a\m) {m\a\n) 



(B7) 



We may now write 



W{n\a\k){k\a\n){n\a\ 



k^n 



^nk 



a\n) ^(n|cr|fc) 



; (en + Cfe)^ 



k^n 



•'nk 



^ ((niacin) - 2{n\a^\n){n\a\n) + {n\a\n)') + ^1 ^"'^'^^ 



2 f— ^ tOnk 



\cr^\n) - {k\a\n){n\a\n)] 



4EE 

k^7i m^n 



{n\a\k) {k\a\m) (m|cr|n) 

^nk^7im 
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With simple algebra we obtain: 



k^n ^"'^ 



k^n 



U>nk 



^^EE 

k^n m^n 
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N 








10 


5.7831903186 


14.682015349 


26.374703996 


20 


5.7831860077 


14.681971199 


26.374617574 


30 


5.7831859657 


14.681970679 


26.374616505 


40 


5.7831859633 


14.681970647 


26.374616438 


50 


5.7831859630 


14.681970643 


26.374616430 


60 


5.7831859630 


14.681970642 


26.374616428 


exact 


5.7831859629 


14.681970642 


26.374616427 



TABLE I: Three lowest energy eigenvalues of a circular mem- 
brane of unit radius obtained using the Collocation Mapping 
Method with different number of elements. The conformal 
map is approximated with the Taylor series up to order z^"^ . 



N 


Eo 


-El, 2 


E-A.i 


20 


5.7833478471 


14.683030231 


26.376321659 


40 


5.7831962133 


14.682036989 


26.374723431 


60 


5.7831879924 


14.681983747 


26.374637569 


80 


5.7831866056 


14.681974789 


26.374623117 


100 


5.7831862262 


14.681972340 


26.374619167 


exact 


5.7831859629 


14.681970642 


26.374616427 



TABLE 11: Three lowest energy eigenvalues of a circular mem- 
brane of unit radius obtained using the Collocation Colloca- 
tion Method with different grid size. 





ny) 


d 


PTo 


PTi 


PT2 


PTa 


exact 


(1 


1) 


1 


4.93480 


5.66472 


5.76740 


5.78118 


5.78319 


(1 


2) 


2 


12.3370 


14.4197 


14.6695 


14.6844 


14.68197 


(2 


1) 




12.3370 


14.4197 


14.6695 


14.6844 


14.68197 


(2 


2) 


1 


19.7392 


25.0228 


26.1714 


26.3593 


26.37462 


(1 


3) 


2 


24.6740 


26.4137 


26.4007 


26.3785 


26.37462 


(3 


1) 




24.6740 


30.2612 


30.7573 


30.5467 


30.47126 


(2 


3) 


2 


32.0762 


40.3174 


41.4936 


41.1156 


40.70647 


(3 


2) 




32.0762 


40.3174 


41.4936 


41.1156 


40.70647 


(1 


4) 


2 


41.9458 


47.7922 


48.5522 


48.9791 


49.21846 


(4 


1) 




41.9458 


47.7922 


48.5522 


48.9791 


49.21846 


(3 


3) 


1 


44.4132 


55.1980 


57.6642 


57.9679 


57.58294 


(2 


4) 


2 


49.3480 


57.1883 


57.7144 


57.6217 


57.58294 


(4 


2) 




49.3480 


66.4508 


71.4408 


72.0826 


70.85000 


(3 


4) 


2 


61.6850 


76.3544 


77.2464 


76.2107 


76.93893 


(4 


3) 




61.6850 


76.3544 


77.2464 


76.2107 


76.93893 


(1 


5) 


2 


64.1524 


72.9892 


71.1717 


70.3712 


70.85000 


(5 


1) 




64.1524 


72.5674 


73.5746 


74.6076 


74.88701 


(5 


2) 


2 


71.5546 


89.4650 


96.8822 


99.3590 


95.27757 


(2 


5) 




71.5546 


89.4650 


96.8822 


99.3590 


95.27757 


(4 


4) 


1 


78.9568 


97.3221 


97.4219 


96.4281 


98.72627 



TABLE III: Energy of the twenty lowest states of the circular 
membrane calculated with perturbation theory to different 
orders. 



26 



n 




a 


1 


" inn 


a - 


1 

- ^nn 




a. = 


1 

- snn 


PT 


1 


-11.9314 


-11 


9828 


-11.9957 


-11 


9989 


-11.9997 


-11 


9999 


-12 


2 


-62.2270 


-62 


9170 


-63.0930 


-63 


1372 


-63.1482 


-63 


1510 


-63.15197799 


3 


-7.74588 


-7 


6738 


-7.6556 


-7 


6510 


-7.6498 


-7 


6495 


-7.649505501 


4 


-12.8722 


-12 


2228 


-12.0560 


-12 


0139 


-12.0034 


-12 


0007 


-12 


5 


-207.2669 


-212 


6755 


-214.1045 


-214 


4669 


-214.5578 


-214 


5805 


-214.5883271 


6 


-13.27723 


-13 


2145 


-13.1989 


-13 


1950 


-13.1940 


-13 


1938 


-13.19388958 


7 


-25.0501 


-24 


8644 


-24.8184 


-24 


8069 


-24.8040 


-24 


8033 


-24.80336862 


8 


3.2132 


8 


2394 


9.5815 


9 


9228 


10.0085 


10 


0299 


10.03674606 


9 


-529.3273 


-558 


6905 


-567.0588 


-569 


2278 


-569.7751 


-569 


9123 


-569.9586259 


10 


-22.0706 


-21 


9810 


-21.9587 


-21 


9531 


-21.9517 


-21 


9513 


-21.95180463 


11 


-12.8590 


-12 


2085 


-12.0512 


-12 


0123 


-12.0026 


-12 


0002 


-12 


12 


-44.8505 


-44 


7072 


-44.6716 


-44 


6627 


-44.6605 


-44 


6599 


-44.66052201 


13 


46.7717 


74 


1811 


82.1277 


84 


1959 


84.7183 


84 


8492 


84.89208802 


14 


-39.7638 


-39 


5706 


-39.5237 


-39 


5120 


-39.5091 


-39 


5084 


-39.50940841 


15 


30.6629 


33 


2746 


33.8759 


34 


0227 


34.0591 


34 


0682 


34.0699668 


16 


-1109.8341 


-1223 


2138 


-1259.7952 


-1269 


6749 


-1272.1960 


-1272 


8294 


-1273.042285 


17 


-33.6433 


-33 


4920 


-33.4543 


-33 


4449 


-33.4425 


-33 


4419 


-33.44310862 


18 


-70.2925 


-70 


1087 


-70.0628 


-70 


0513 


-70.0484 


-70 


0477 


-70.04918774 


19 


102.5451 


207 


3338 


242.1806 


251 


6572 


254.0794 


254 


6883 


254.8899234 


20 


-12.5604 


-12 


1285 


-12.0299 


-12 


0059 


-11.9999 


-11 


9984 


-12 


21 


-70.6796 


-70 


5877 


-70.5651 


-70 


5595 


-70.5581 


-70 


5577 


-70.55997503 


22 


124.8295 


134 


9087 


137.0858 


137 


6039 


137.7317 


137 


7635 


137.7717825 


23 


-2007.8189 


-2335 


7052 


-2460.1854 


-2496 


0792 


-2505.4160 


-2507 


7742 


-2508.565211 


24 


-47.9136 


-47 


6534 


-47.5887 


-47 


5725 


-47.5685 


-47 


5675 


-47.56994783 


25 


-101.2945 


-101 


0139 


-100.9438 


-100 


9262 


-100.9219 


-100 


9208 


-100.9237065 


26 


115.2123 


410 


2533 


528.4606 


562 


9517 


571.9490 


574 


2230 


574.980025 


27 


-53.7499 


-53 


6334 


-53.6067 


-53 


6002 


-53.5986 


-53 


5981 


-53.60145365 


28 


50.1431 


51 


1784 


51.3848 


51 


4331 


51.4450 


51 


4480 


51.44555258 


29 


-106.3989 


-106 


3559 


-106.3450 


-106 


3423 


-106.3416 


-106 


3414 


-106.3455495 


30 


286.6324 


322 


1690 


329.4268 


331 


0946 


331.5016 


331 


6027 


331.632198 


31 


-3256.1580 


-4011 


4680 


-4354.9098 


-4463 


3391 


-4492.3946 


-4499 


7924 


-4502.275663 


32 


-64.8743 


-64 


4402 


-64.3324 


-64 


3055 


-64.2988 


-64 


2971 


-64.30174616 


33 


-12.0909 


-12 


0078 


-11.9974 


-11 


9955 


-11.9950 


-11 


9949 


-12 


34 


-96.1084 


-96 


0849 


-96.0799 


-96 


0787 


-96.0784 


-96 


0783 


-96.08384485 


35 


167.0112 


169 


8399 


170.2940 


170 


3920 


170.4156 


170 


4214 


170.4178176 


36 


-137.8993 


-137 


4572 


-137.3467 


-137 


3190 


-137.3121 


-137 


3104 


-137.3156766 


37 


-0.2977 


647 


4859 


970.7556 


1074 


8721 


1102.9016 


1110 


0463 


1112.43464 


38 


-147.5040 


-147 


5073 


-147.5074 


-147 


5074 


-147.5073 


-147 


5073 


-147.5143637 


39 


516.9942 


627 


1522 


649.2405 


654 


1087 


655.2782 


655 


5675 


655.6568115 


40 


-67.3226 


-67 


3763 


-67.3928 


-67 


3972 


-67.3982 


-67 


3985 


-67.40630759 


41 


67.0772 


66 


9190 


66.8387 


66 


8163 


66.8105 


66 


8091 


66.8009496 


42 


-4871.0690 


-6337 


4062 


-7134.0401 


-7415 


2895 


-7493.8895 


-7514 


1430 


-7520.958743 


43 


-141.9693 


-142 


0124 


-142.0235 


-142 


0262 


-142.0270 


-142 


0271 


-142.0360122 


44 


-84.5423 


-83 


8471 


-83.6747 


-83 


6317 


-83.6210 


-83 


6183 


-83.62637707 


45 


359.2556 


367 


9325 


369.0389 


369 


2503 


369.2992 


369 


3111 


369.3064766 


46 


-330.0113 


-179 


4676 


-179.2960 


-179 


2531 


-179.2424 


-179 


2397 


-179.2486248 


47 


-180.1538 


837 


1280 


1575.4900 


1844 


8179 


1920.6443 


1940 


2181 


1946.790496 


48 


-11.4814 


-11 


8521 


-11.9544 


-11 


9805 


-11.9871 


-11 


9887 


-12 


49 


-194.2969 


-194 


3597 


-194.3737 


-194 


3770 


-194.3779 


-194 


3781 


-194.3894759 


50 


785.2942 


1078 


8172 


1140.1836 


1153 


1570 


1156.2105 


1156 


9613 


1157.199583 



TABLE IV: Corrections to the energy of a particle in a box of square 2 for deformations of the square, {En ~ E'^^)/o? , for 
different values of a for the first 50 energy levels. In the last column are the results of Perturbation Theory. 
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TABLE V: First order perturbative corrections (divided by the corresponding energy in the square) for the first twenty states 
of a square drum subject to an arbirtrary (weak) deformation. For the states marked with f the degeneracy is hfted by the 
perturbation. 
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A = 1/100 


A = 1/20 


(fc, n) 


PTo 


PT2 


CCM 


PTo 




CCM 


(0,1) 


5.78434 


5.78319 


5.78319 


5.81210 


5.78304 


5.78325 


(1,1)^ 


14.6849 


14.6805 


14.6805 


14.7554 


14.6447 


14.646 


(1,1)^ 


14.6849 


14.6834 


14.6834 


14.7554 


14.7185 


14.7183 


(2,1) 


26.3799 


26.3746 


26.3746 


26.5065 


26.3740 


26.3734 


(2,1) 


26.3799 


26.3746 


26.3746 


26.5065 


26.3740 


26.3741 


(0,2) 


30.4774 


30.4713 


30.4713 


30.6236 


30.4705 


30.4739 


(3,1) 


40.7146 


40.7065 


40.7065 


40.9100 


40.7054 


40.7056 


(3,1) 


40.7146 


40.7065 


40.7065 


40.9100 


40.7054 


40.7056 


(1,2) 


49.2283 


49.2135 


49.2136 


49.4645 


49.0936 


49.0991 


(1,2) 


49.2283 


49.2234 


49.2234 


49.4645 


49.3409 


49.3415 


(4,1) 


57.5945 


57.5829 


57.5830 


57.8709 


57.5815 


57.582 


(4,1) 


57.5945 


57.5829 


57.5830 


57.8709 


57.5815 


57.582 


(2,2) 


70.8642 


70.8500 


70.8500 


71.2042 


70.8482 


70.8389 


(2,2) 


70.8642 


70.8500 


70.8500 


71.2042 


70.8482 


70.8501 


(0,3) 


74.9020 


74.8870 


74.8871 


75.2614 


74.8851 


74.9035 


(5,1) 


76.9543 


76.9389 


76.9390 


77.3236 


76.9370 


76.9378 


(5,1) 


76.9543 


76.9389 


76.9390 


77.3236 


76.9370 


76.9378 


(3,2) 


95.2966 


95.2776 


95.2776 


95.7540 


95.2752 


95.2734 


(3,2) 


95.2966 


95.2776 


95.2776 


95.7540 


95.2752 


95.2736 


(6,1) 


98.7460 


98.7263 


98.7263 


99.2199 


98.7238 


98.7251 


(6,1) 


98.7460 


98.7263 


98.7263 


99.2199 


98.7238 


98.7251 


(l,3)t 


103.520 


103.489 


103.489 


104.017 


103.237 


103.253 


(l,3)t 


103.520 


103.510 


103.510 


104.017 


103.757 


103.763 


(4,2) 


122.452 


122.428 


122.428 


123.040 


122.425 


122.422 


(4,2) 


122.452 


122.428 


122.428 


123.040 


122.425 


122.422 


(7,1) 


122.932 


122.908 


122.908 


123.522 


122.904 


122.908 


(7,1) 


122.932 


122.908 


122.908 


123.522 


122.904 


122.908 


(2,3) 


135.048 


135.021 


135.021 


135.696 


135.017 


134.978 


(2,3) 


135.048 


135.021 


135.021 


135.696 


135.017 


135.025 


(0,4) 


139.068 


139.040 


139.041 


139.735 


139.037 


139.098 


(8,1) 


149.483 


149.453 


149.453 


150.200 


149.449 


149.451 


(8,1) 


149.483 


149.453 


149.453 


150.200 


149.449 


149.451 


(5,2) 


152.272 


152.241 


152.241 


153.002 


152.237 


152.238 


(5,2) 


152.272 


152.241 


152.241 


153.002 


152.237 


152.238 


(3,3) 


169.429 


169.395 


169.396 


170.242 


169.391 


169.383 


(3,3) 


169.429 


169.395 


169.396 


170.242 


169.391 


169.384 


(l,4)t 


177.556 


177.503 


177.503 


178.408 


177.070 


177.108 


(l,4)t 


177.556 


177.539 


177.539 


178.408 


177.962 


177.982 


(9,1) 


178.373 


178.337 


178.338 


179.229 


178.333 


178.335 


(9,1) 


178.373 


178.337 


178.338 


179.229 


178.333 


178.335 



TABLE VI: First 40 eigenvalues of the Robnik billiard with A = 1/100 and A = 1/20. The CCM results are obtained with a 
grid with A'^ = 100. The symbol f is used to highlight the states where first order perturbation theory predicts a lifting of the 
degeneracy. 
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octagon 


nonagon 


decagon 


{k, 


n) 


PTo 


PTi 


CCM 


PTo 


PTi 


CCM 


PTo 


PTi 


CCM 


(0 


1) 


6.48669 


6.48505 


6.48493 


6.32407 


6.32314 


6.32309 


6.21258 


6.21202 


6.21200 


(1 


1) 


16.468 


16.4581 


16.4561 


16.0551 


16.0495 


16.0486 


15.7721 


15.7687 


15.7682 


(1 


1) 


16.468 


16.4581 


16.4561 


16.0551 


16.0495 


16.0486 


15.7721 


15.7687 


15.7682 


(2 


1) 


29.583 


29.5530 


29.5406 


28.8413 


28.8241 


28.8185 


28.3329 


28.3224 


28.3197 


(2 


1) 


29.583 


29.5530 


29.5406 


28.8413 


28.8241 


28.8186 


28.3329 


28.3224 


28.3197 


(0 


2) 


34.178 


34.1407 


34.1245 


33.3211 


33.2994 


33.2920 


32.7337 


32.7204 


32.7166 


(3 


1) 


45.6583 


45.5912 


45.5298 


44.5136 


44.4747 


44.4501 


43.7289 


43.7050 


43.6936 


(3 


1) 


45.6583 


45.5912 


45.5298 


44.5136 


44.4747 


44.4501 


43.7289 


43.7050 


43.6937 


(1 


2) 


55.2057 


55.1197 


55.0498 


53.8217 


53.7708 


53.7391 


52.8729 


52.8412 


52.8254 


(1 


2) 


55.2057 


55.1197 


55.0498 


53.8217 


53.7708 


53.7391 


52.8729 


52.8412 


52.8255 


(4 


1) 


64.5877 


62.6348 


62.5959^ 


62.9685 


62.8946 


62.7662 


61.8584 


61.8128 


61.7682 


(4 


1) 


64.5877 


66.2878 


66.2775+ 


62.9685 


62.8946 


62.7662 


61.8584 


61.8128 


61.7682 


(2 


2) 


79.4687 


79.3097 


79.0194 


77.4763 


77.3810 


77.2694 


76.1106 


76.0505 


75.9987 


(2 


2) 


79.4687 


79.3097 


79.0199 


77.4763 


77.3810 


77.2696 


76.1106 


76.0505 


75.9988 


(0 


3) 


83.9968 


83.8282 


83.5906 


81.8909 


81.7892 


81.6829 


80.4473 


80.3828 


80.3302 


(5 


1) 


86.2983 


86.0854 


85.8297 


84.1347 


84.0094 


84.0635 


82.6516 


81.0526 


81.0237+ 


(5 


1) 


86.2983 


86.0854 


85.8297 


84.1347 


84.0094 


84.0635 


82.6516 


84.0950 


84.0826+ 


(3 


2) 


106.868 


106.609 


106.713 


104.189 


104.032 


102.757 


102.352 


102.252 


102.055 


(3 


2) 


106.868 


106.609 


106.713 


104.189 


104.032 


102.757 


102.352 


102.252 


102.055 


(6 


1) 


110.736 


110.405 


110.452 


107.960 


107.763 


108.953 


106.057 


105.934 


105.808 


(6 


1) 


110.736 


110.405 


110.452 


107.960 


107.763 


108.953 


106.057 


105.934 


105.808 


(1 


3) 


116.090 


115.812 


114.883 


113.179 


113.009 


112.684 


111.184 


111.075 


110.926 


(1 


3) 


116.090 


115.812 


114.883 


113.179 


113.009 


112.684 


111.184 


111.075 


110.927 


(4 


2) 


137.321 


133.003 


132.733+ 


133.878 


133.641 


133.313 


131.518 


131.366 


131.401 


(4 


2) 


137.321 


137.373 


137.954+ 


133.878 


133.641 


133.313 


131.518 


131.366 


131.401 


(7 


1) 


137.859 


137.373 


137.954 


134.403 


134.113 


133.485 


132.033 


131.851 


131.925 


(7 


1) 


137.859 


140.865 


140.802 


134.403 


134.113 


133.485 


132.033 


131.851 


131.926 


(2 


3) 


151.446 


151.031 


150.459 


147.649 


147.391 


147.814 


145.046 


144.879 


144.205 


(2 


3) 


151.446 


151.031 


150.461 


147.649 


147.391 


147.815 


145.046 


144.879 


144.205 


(0 


4) 


155.954 


155.530 


151.944 


152.044 


151.780 


150.877 


149.364 


149.192 


148.819 


(8 


1) 


167.633 


165.334 


168.290+ 


163.431 


163.022 


162.822 


160.550 


160.292 


160.168 


(8 


1) 


167.633 


168.569 


168.776+ 


163.431 


163.022 


162.822 


160.550 


160.292 


160.817 


(5 


2) 


170.761 


170.216 


168.776 


166.480 


166.143 


166.165 


163.545 


160.293 


160.817+ 


(5 


2) 


170.761 


170.216 


169.209 


166.480 


166.143 


166.165 


163.545 


166.360 


166.328+ 


(3 


3) 


190.002 


189.422 


189.910 


185.238 


184.875 


183.416 


181.973 


181.735 


181.436 


(3 


3) 


190.002 


189.422 


189.910 


185.238 


184.875 


183.416 


181.973 


181.735 


181.437 


(1 


4) 


199.116 


198.519 


192.516 


194.124 


193.008 


193.117 


190.702 


190.454 


186.961 


(1 


4) 


199.116 


198.519 


192.516 


194.124 


193.748 


193.117 


190.702 


190.454 


186.961 


(9 


1) 


200.032 


199.109 


204.554 


195.017 


193.748 


193.659+ 


191.579 


191.226 


194.384 


(9 


1) 


200.032 


199.109 


204.554 


195.017 


195.913 


195.691+ 


191.579 


191.226 


194.384 



TABLE VII: First 40 eigenvalues of the octagon, nonagon and decagon billiards. The symbol f is used to highlight the states 
where first order perturbation theory predicts a lifting of the degeneracy. The third column for each polygon displays the 
numerical result obtained using the CCM with a grid corresponding to A*' = 100. 



